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Limiters are nonlinear hybridization techniques that are used to preserve positivity and 
monotonicity when numerically solving hyperbolic conservation laws. Unfortunately, the 
original methods suffer from the truncation-error being I s * order accurate at all extrema 
despite the accuracy of the higher-order method jU [21 [31 H] . To remedy this problem, higher- 
order extensions were proposed that relied on elaborate analytic and geometric constructions 
[S, El [71 E] • Since extremum-preserving limiters are applied only at extrema, additional com- 
putational cost is negligible. Therefore, extremum-preserving limiters ensure higher-order 
spatial accuracy while maintaining simplicity. This report presents higher-order limiting for 
(z) computing van Leer slopes and (ii) adjusting parabolic profiles. This limiting preserves 
monotonicity and accuracy at smooth extrema, maintains stability in the presence of discon- 
tinuities and under-resolved gradients, and is based on constraining the interpolated values 
at extrema (and only at extrema) by using nonlinear combinations of 2 nd derivatives. The 
van Leer limiting can be done separately and implemented in MUSCL (Monotone Upstream- 
centered Schemes for Conservation Laws) [2j or done in concert with the parabolic profile 
limiting and implemented in PPM (Piecewise Parabolic Method) [91 [TOj . The extremum- 
preserving limiters elegantly fit into any algorithm which uses conventional limiting tech- 
niques. Limiters are outlined for scalar advection and nonlinear systems of conservation 
laws. This report also discusses the 4 th order correction to the point-valued, cell-centered 
initial conditions that is necessary for implementing higher-order limiting. The material 
herein complements Colella and Sekora [TT]. Lastly, there is no guarantee that extremum- 
preserving limiters preserve positivity. To ensure this property, one should combine the 
limiting with FCT (Flux-Corrected Transport) [3]. 
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1 Algorithms for Scalar Advection 

Consider the following scalar equation in one spatial dimension: 

da da . . 

— — V A— — = 0, A = constant. (1) 
at ox 

At time-step n, the average-valued, cell-centered quantity a over a finite volume of length 
h = Ax is: 

1 r(i+h)h 

a i ~ r / a(aj,nAt)cfe. (2) 
^ V(i-|)h 

MUSCL/PPM are conservative finite volume methods that are used to compute a" +1 : 

XAt 



a? +1 = a? 



Ax 



(v? - <-1) • (3) 



where is the average of a linear/parabolic interpolant over the interval swept out by 
the characteristics crossing the cell face at {i + \)h and is given by: 



where a = G [0, 1] is the CFL number and aj(x) is the linear/parabolic interpolant, 
such that x E [(i — \)h, (i + \)h\. 

There are three variations of Godunov-type methods for which limiters can be implemented: 

1. van Leer Limiter in MUSCL 

2. van Leer Limiter + Parabolic Profile Limiter in PPM 

3. Parabolic Profile Limiter in PPM 

Each of these algorithm variations are discussed below. 

1.1 MUSCL 

1. van Leer limit the differences (Acij), giving 2 nd order results. Apply the corresponding 
boundary conditions. 

2. Use the van Leer limited differences to compute 4 th order differences: 

A 4 ai = - |(a j+ i - |Aa !+! ) - (a,_i + ^Aa^)) . (5) 
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3. Employ piecewise linear reconstruction by computing spatially extrapolated face-centered 
values at the low and high (left and right) edges of cells: 

1 

Oi,± = aj + - (±1 - a) A 4 aj. (6) 

1.2 PPM 

1. van Leer limit the differences (Aa*), giving 2 nd order results. Apply the corresponding 
boundary conditions. 

2. Employ either 4 th or 6* h order piecewise parabolic reconstruction: 

• Method 1: use the van Leer limited differences and compute spatially extrapolated 
face-centered values at the low and high (left and right) edges of cells: 



4 

ai,+ 


= ^ (ai+i + ai) - ^ (Aa i+ i - Acii) , 






(7) 


4 


= ^ (^ + ai_i) - ^ (Aoj - Aai_i) , 






(8) 


v+ 


= 4 a ii+ - ^ (3(Aoi + i - Aai) - (Aai +2 


- Aai_ 


-i)), 


(9) 


V 


= \ - - (3(Aai - Aaj-i) - (Aai+i 


- Aai_ 


-2)), 


(10) 




= a± - ai. 






(11) 



• Method 2: employ either 4 th or 6 th order piecewise parabolic reconstruction with- 
out using the van Leer limited differences in Step 1: 



4 a i+ i = — (oj+i + a t ) - (a i+2 + Oi-i) , (12) 
37 , 8 , ,1 



7 1 

**+3 ~~ 12 + ai > ~ 12 

6 Oi + i = ^ (fli+i + fli) - ^ (a i+2 + Oi_i) + ^ (oi +3 + ai_ 2 ) , (13) 
a± = a i± i - ai. (14) 

It is important to note that when Acti = A c aj (centered difference) the two for- 
mulations for the piecewise parabolic reconstruction are identical. 

3. Limit the parabolic profile (a±). 

4. Use the PPM predictor values to reconstruct the parabolic profile: 

a± = ai + a± + ^(±(a_ — a+)) — (a- + «+)(3 — 2a). (15) 
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1.3 Update Solution for MUSCL or PPM 



1. Compute fluxes: 




Aa+ A > 0, 
Aa_ A < 0. 



(16) 



2. Use the divergence of the fluxes to update the solution: 



n+l 



a 




Ax V l+i 3 




(17) 



- - 



A > 



A < 0. 



(18) 



1.4 Conventional Limiters 

To complete the specification of the MUSCL and PPM schemes, one defines the conventional 
limiters used to constrain the interpolated profiles within each cell. 

1.4.1 Conventional van Leer Limiter 

Given a sequence of average- valued, cell-centered quantities, the conventional van Leer limiter 
proceeds with the following steps [21 HE HH] : 

1. Compute one-sided and centered differences: 



A c di 

A + CLi 



A-di 



-{a-i+i — cii-i), 



(19) 
(20) 
(21) 



2. Apply the conventional van Leer limiter: 





(23) 



(24) 



One significant defect of this method is the clipping of extremum when A_<2j A + aj < 0. 
This clipping sets Aaj — * as a precautionary measure for suppressing spurious oscillations. 
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How one arrives at the formula for the conventional van Leer limiter can be understood by 
considering the following example. Assume that one is not at an extremum, ^\ih < 0, and 
|A_dj| > |A + Oj|. The value of a((i + \)h) at face-center (i + is approximated by: 

a{{i + \)h) = a((i + l)h) - + •■■ 

a((i + l)h) = a{ih) + ffU + ..• 

Clearly, a((i + \)h) > a((i + l)h) and a((i + \)h) « a(ih) + \%\ih- Therefore: 

h^\ ih >2(a((i + l)h)-a(ih)), (26) 
and one arrives at Eq. [221 However, at extrema ^ — > 0, which leads one to Eq. [2H 

1.4.2 Conventional Parabolic Profile Limiter 

Given the higher-order reconstruction of a± such that a± = a± — a,i, the conventional 
parabolic profile limiter proceeds with the following steps 0, [10] : 

1. Adjust a± according to the following cases: 

a± — > > 0, 

a + -> -2a_ > 4a 2 _, (27) 

a_ — >■ — 2a + a 2 _ > Aa 2 ^. 

2. Reconstruct a± given the adjusted values for a±. 

One significant defect of this method is the adjustments to a± make the reconstruction a 
monotone profile as a precautionary measure for suppressing spurious oscillations. However, 
this constraint is more restrictive than is required to preserve monotonicity [T3] . 

The adjustments to a± are derived from the interpolation polynomial described in the orig- 
inal Piecewise Parabolic Method: 

Oj = a_ + a(8±a + a 6 (l — a)), (28) 
5±a = a + — a_ = a + — a_, (29) 

a 6 = 6 ^ai - ^(a + + a_)^j = — 3(a + + a_). (30) 

where a G [0, 1] is the CFL number, which also corresponds to a dimensionless length scale. 
This scale is associated with each grid cell such that the left side of a cell is designated a = 
and the right side of a cell is designated a = 1. Differentiating Oj with respect to a gives: 

— — = 5±a + a 6 (l — 2cr) = (a + — a;_) — (a + + a_)(l — 2cr). (31) 
ad 
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By evaluating da^jda at the left and right sides of the cell: 

= 5±a + a 6 = -2(a+ + 2a_), (32) 

cr^O 

= S±a — ciq = 2(2a + + a_). (33) 



dcii 
da 



<T-t-l 



Maximize with respect to a by setting the derivatives equal to zero and solving the 
resulting equations. One arrives at the following result: 

a + -> -2a_ a -> 0, 
ct_ — > — 2a + <j — > 1. 

2 Extremum-Preserving Limiters 

For smooth solutions away from extrema, the MUSCL scheme is 2 nd order accurate for linear 
advection whereas the PPM is 3 rd order accurate for linear advection and 4 th order accurate 
in the limit of vanishing CFL number. However, the monotonicity constraints at extrema 
reduce the truncation error to 0(h) even at smooth extrema. This reduction in the overall 
accuracy of the method also introduces a non-smooth component to the error. To eliminate 
this problem, one constructs a new limiting scheme at extrema. 

The defect of the standard approach to limiting is most easily seen in the MUSCL limiter. 
Away from extrema, the magnitude of the slope is computed as the minimum of three un- 
divided differences: the centered difference and twice the one-sided differences. In smooth 
regions away from extrema, the centered difference and the one-sided differences all ap- 
proximate h^ such that the minimum is always defined by the centered difference. At 
discontinuities, one of the one-sided differences is typically much smaller than the other two 
differences. Therefore, this one-sided difference is chosen because it leads to a reduction in 
the slope and suppresses oscillations. However, the idea behind this method fails at extrema 
because the derivative vanishes. Furthermore, the one-sided differences have opposite signs 
and this non-constant multiple bounds the centered difference. In the original van Leer and 
PPM limiters, the solution is to simply drop the order of the method to 1 st order. 

In the approach used in this report, one changes the limiters at extrema, and only at extrema, 
by using comparisons of different estimates of the 2 nd derivatives as a basis for whether to 
limit the underlying linear scheme. If the solution is smooth at the extremum, then all 
of the estimates of the 2 nd derivative are comparable and the limiter leaves the underlying 
linear scheme unchanged. Discontinuities, underresolved gradients, and high-wavenumber 
oscillations are detected either by one of the estimates of the 2 nd derivative being much 
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smaller than the others or by the various estimates of the 2 nd derivatives changing sign. 
Either effect triggers a nontrivial limiting of the interpolating function in the cell and a 
resulting suppression of oscillations. 

2.1 Extremum-Preserving van Leer Limiter 

The extremum-preserving van Leer limiter parallels the conventional van Leer limiter: 

1. Compute one-sided and centered differences as well as one-sided differences that are 
an additional spatial step-size away: 

A — Oj = cij_i - aj_ 2 , (35) 

A_aj = aj - aj_i, (36) 

A c aj = ^(oj+i -Oi_i), (37) 

A + ai = a i+ i - aj, (38) 

A ++ a; = a i+2 -a i+ i. (39) 

2. Test for extrema. An extremum is defined if the following condition is satisfied: 

min(A_ a^A+aj, A a* A ++ a;) < 0. (40) 

3. If the above extremum condition is not satisfied, then limiting follows the conventional 
van Leer method: 

Ai im ai = 2min(|A_a;|, |A + a;|), (41) 

S = sign(A c a,), (42) 

Aaj = min(|A c aj|, |Ai im ai|)5. (43) 

4. If the above extremum condition is satisfied, then: 

• Compute one-sided and centered 2 nd derivatives: 

= ^(a(ih)-2a((i-l)h) + a((i-2)h)), (44) 

V 2 c a t = ^ (a((i + l)h) - 2a(ih) + a((i - l)h)) , (45) 

V\a { = ^(a{(i + 2)h)-2a{(i + l)h) + a(ih)). (46) 

• Find the minimum 2 nd derivative over the five-cells in question: 

S 2 = sign(P c 2 a,), (47) 

Vf im ai = mm(\V 2 c a i \,max(S 2 V^a i ,0),m&x(S 2 Vla i ,0)). (48) 
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• Apply the modified van Leer limiter at the extremum: 

f minfCy L ^P 1 2 im «i,2|A_a i |) 5 2 A c a, < 

min(^Cy L ^PL«, ? 2|A + a,|J else 

5 = sign(A c aj), (50) 

Aflj = min(|A c aj|, |A Um aj|)<S. (51) 

V designates a derivative while A designates a difference. Derivatives and differences are 
similar operators that can be transformed back-and-forth when one considers the relevant 
Taylor expansion and multiplies/divides each operator by factors of h. Cvl is a constant 
that is independent of the mesh spacing and as Cvl — > 0, the extremum-preserving van Leer 
limiter reduces to the conventional van Leer limiter. For most calculations, Cvl = 1-25. 

How one arrives at the tighter bound of ^§-^ij m ai for peak height in the extremum-preserving 
van Leer limiter can be understood by considering the following example. Assume that one 
is at a local maximum such that ^\ih < and < 0. The van Leer limiting condition 
bounds the derivative on the high (right) edge of the cell: 

Aa i >2(a i+1 -a i ). (52) 

When looking for a bound on Acij from the low (left) edge of the cell, an extremum near % 
implies that: 

< a((i - l)h) - a((i - 2)h) = ^l M)fc + °(^)- ( 53 ) 
Therefore, one arrives at the following bound on Aaj: 

a i da. , da, 3h 2 d 2 a 3h 2 d 2 a 

Aa, = h-U + 0(*«) = h-\ { ,_, )h + + 0(fc» ) > + 0(*«). (54) 

These bounds are summarized in the following inequality: 

0>Aa, t >C VL — — , C yi >l, (55) 

where the 3/2 comes from the condition that (i — |)/i is the nearest face-centered point to 
the cell-centered point ih for which one can unequivocally assert that h^\^_^ h > 0(h). 



2.2 Extremum-Preserving Parabolic Profile Limiter 

If the piecewise parabolic reconstruction is done without using van Leer limited differences, 
then one has to include an additional step that limits extrema at cell faces such that a i+ i 
lies between adjacent cell averages. Van Leer limiting automatically enforces this constraint. 
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• Test for extremum at cell faces: 

(a i+ i - ai){a i+1 - a i+ i) < 0. (56) 



If no extremum is found, then proceed without any adjustment to a 



l . 



• If an extremum is found, compute one-sided and centered 2 nd derivatives: 

V 2 _a l+ i = ^ (a((i + l)h) - 2a(ih) + a((i - l)h)) , (57) 
2 h 



V 2 c a i+h = - [a((i + l)h) - 2a((i + -)h) + a(ih) j , (58) 

V 2 + a l+ i = —(a((i + 2)h)-2a((i + l)h)+a(ih)). (59) 

• Find the minimum difference with respect to o i+ i that is greater than zero: 

Sf H = sign^Vi), (60) 
A 2 im«i+i = ^x(mm(Cp PM S 2 + iV 2 _a i+ i,Sf + iV 2 c a i+ i, 

C PPM S 2 H Vla i+h ),0). (61) 

• Adjust a i+ i according to: 

a i+ i = ^(a i+ i + ai) - ^Vl m a i+ i. (62) 
The above formulas are arrived at by considering the centered 2 nd derivative: 

V " a ^ = (hj2j (° (( * + l)k) ~ 2a(( ' + \ )k) + a(ik) ) ■ (63) 

The spatial resolution oih/2 enters because one is considering differences between a(ih), a((i+ 
\)h), a((i + l)h). By substituting in the following expressions: 

a, = a(ih) -^D 2 c ai + 0(h 4 ) : (64) 

a i+1 = a((i + l)h)-^V 2 c a t+1 + 0(h A ), (65) 
which convert point- valued to average- valued quantities, one arrives at: 

^V 2 c a l+k = a i+1 - 2a((i + l -)h) + a t - + 0{h A ). (66) 



Rearranging the above equation gives the expression: 

3 



V 2 a l+ , = — ( a i+1 - 2a((i + -)h) + a t ) . (67) 
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After finding the minimum 2 nd derivative V 2 m a i+ i , one can again use the above equation to 
assign a value to a i+ i: 

a i+\ = \ ( a i+i + a i) ~ \ v lm a i+\- ( 68 ) 
Lastly, a± = a i± i and one proceeds with limiting a±. 

Now given the higher-order reconstruction of a± such that a± = a± — a^ adjust a± according 
to the following cases: 

1. a+a- > ] (a i+1 - ai)(a,i - di-ii <0 

• Compute the 2 nd derivative with respect to o from the the interpolation polyno- 
mial described in the original Piecewise Parabolic Method: 

2 d 2 di , 

T> PPM ai = —— = -2a 6 = Q(a + + a^). (69) 
era 

It is important to note that this 2 nd derivative with respect to a is also a 2 nd order 
difference for a,. Furthermore, this difference represents the maximum difference 
that can occur across a cell given a parabolic profile. 

• Compute one-sided and centered 2 nd order derivatives: 

V 2 _di = (a(ih) - 2a((i - l)h) + a((i - 2)h)) , (70) 
h 2 

V 2 c a t = ^(a((i + l)h)-2a(ih)+a((i-l)h)), (71) 

VUi = ^r(a((i + 2)h)-2a((i + l)h) + a(ih)). (72) 

• Find the minimum difference that is greater than zero, given the 2 nd derivatives 
over the five-cells in question as well as the above difference that was derived from 
the interpolation polynomial: 

Sppm = sign(£>f> PM a;), (73) 
"^fiiiA = max ( m i n (<SppM'E > PPM a ii CppMSp PM V 2 _ai^ 

CppM<Sp PM T>' 2 ai, CppMSppj^V^cii), 0). (74) 

• Adjust a± according to T> 2 m cii and V 2 PPM ai\ 

I ^PPM a i\ 

2. a% > 4a 2 _ 
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• Compute the maximum value of a_ over a given cell: 

-a 2 

= -. (76) 

This formula is arrived at by averaging the interpolation polynomial from the 
original Piecewise Parabolic Method on the left side of cell v. 

f- = - r<m (77) 

a Jo 

= a_ + |(M+(l-f aejj (78) 
= a_ + -(a + -a_) - y I 1 - y 1 (a + + a_). (79) 
Maximize /_ with respect to er: 

= °— = 2KTZ)' (80) 



Therefore: 



a max = /-(a max ) - a, = ^ + (81) 

• Preserve monotonicity by ensuring that the following condition is just satisfied: 

I- < di-i -Oi=>I- = Oi_i -Oi = a max . (82) 

• Solve the above equation for a + : 

S~ = sign(a;_), (83) 
a+ = -2/_ -2S- (I 2 -a_/_) 1/2 . (84) 



3. o?L > <±ol\ 



Compute the maximum value of a + over a given cell: 

-a 2 

a T X = 77 — -v ( 85 ) 

+ 4(a+ + a„) 

This formula is arrived at by averaging the interpolation polynomial from the 
original Piecewise Parabolic Method on the right side of cell v. 

U = ~ f <0i (86) 

°" Jl-a 

= cl + - °- (s±a -(l- ^-a 6 ) ) (87) 



a . . 3cr / ' 2a\ 

= a+- - («+-«_)- — II- — 1 («+ + «_). 



(88) 
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Maximize /+ with respect to a: 



df+ 2a + + a_ , , 

1 = °^- = ^)' (89) 



Therefore: 



-a 2 



< ax = /+(<w) -a t = (90) 



• Preserve monotonicity by ensuring that the following condition is just satisfied: 

1+ > Oi+i -cn=^I+ = a l+1 -ai = a+ ax . (91) 

• Solve the above equation for a_: 

S + = sign(a + ), (92) 
a_ = -21+ -2S + (I 2 -«+/+) 1/2 • (93) 

4. Reconstruct a± given the adjusted values for a±. 



3 Numerical Results 

Results are presented for ID scalar advection. The standard test problems were employed 
to demonstrate improvements in accuracy [H]. The following parameters were used in all 
test problems: 

a = 0.2, are [0,1], t = 10, A = 1. (94) 
Periodic boundary conditions were used with the following number of ghost cells: 

A th Order Reconstruction =^ 3 ghost cells per side (95) 
6 th Order Reconstruction =^ 4 ghost cells per side. (96) 

The test problems were defined by the following initial conditions: 

Gaussian Wave: a = exp (— 256(x — 0.5) 2 ) (97) 

Senu-Circle Wave: a = < <°' 252 " <* " °^ ^ < * < °' 75 (98) 

^ else V ; 

. 1 0.25<x<0.75 . , 

Square Wave: a = < (99) 

else 

The conventional initialization for 2 nd order accurate methods is obtained by approximating 
the average over a cell by the value at the center of the cell, i.e., the midpoint rule for 
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integrals. However, that initialization is only 2 order accurate, and less accurate than 
what one expects for PPM when applied to the advection equation. For that reason, one 
uses a 4 th order accurate approximation to the cell average, following [12] : 



This 4 order correction is only administered when defining the initial conditions. The 
corresponding boundary conditions are applied and one begins advancing the temporal loop. 
For smooth solutions away from extrema, PPM is 3 rd order accurate for linear advection and 
4 th order accurate in the limit of vanishing CFL number o —>■ 0. The following definitions 
for the n- norms and convergence rates are used throughout this note. Given the numerical 
solution dF at resolution r and the analytic solution u, the error at a given point i is: 



cii = a(ih) + — V^cti, V\a, 



- {a{i + l)h) - 2a(ih) + a{{i - l)h)) . (100) 



V V 

e t = a t -Ui. 



(101) 



The 1-norm and oo-norm of the error are: 




L 



oo 




(102) 



The convergence rate is measured using Richardson extrapolation: 



R. 



In (L n (e r )/L n (e r+1 )) 
In (Ax r /Ax r+1 ) 



(103) 



14 



M Sekora, P Colella 



3.1 4 th Order Reconstruction 



A th Order, Gaussian Wave 




N ce ll 


u 


Ri 


Loo 


Roo 


No Limiter 


32 


8.0E-2 




3.5E-1 






64 


2.6E-2 


1.6 


1.5E-1 


1.2 




128 


3.2E-3 


3.0 


2.8E-2 


2.5 




256 


3.1E-4 


3.4 


3.0E-3 


3.2 


PPM Limiting C PPM = 


32 


7.5E-2 




4.7E-1 






64 


2.7E-2 


1.5 


2.6E-1 


0.9 




128 


7.8E-3 


1.8 


9.9E-2 


1.4 




256 


1.3E-3 


2.6 


3.1E-2 


1.7 


PPM Limiting C PPM = 1-25 


32 


5.5E-2 




3.7E-1 






64 


1.6E-2 


1.8 


1.4E-1 


1.4 




128 


3.2E-3 


2.3 


2.8E-2 


2.4 




256 


3.1E-4 


3.4 


3.0E-3 


3.2 


VL+PPM Limiting C VL , C PPM = 1-25 


32 


5.7E-2 




3.8E-1 






64 


1.6E-2 


1.8 


1.5E-1 


1.4 




128 


3.3E-3 


2.3 


2.8E-2 


2.4 




256 


3.1E-3 


3.4 


3.0E-3 


3.2 
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4 th Order, Semi-Circle Wave 




N ceU 


L x 


Ri 




Roo 


No Limiter 


32 


1.2E-2 




5.2E-2 






64 


5.9E-3 


1.0 


4.1E-2 


0.3 




128 


2.6E-3 


1.2 


3.2E-2 


0.3 




256 


1.1E-3 


1.2 


2.5E-2 


0.3 


PPM Limiting C PPM = 


32 


8.1E-3 




4.1E-2 






64 


4.7E-3 


0.8 


3.4E-2 


0.3 




128 


2.1E-3 


1.2 


2.7E-2 


0.3 




256 


9.0E-4 


1.2 


2.1E-2 


0.3 


PPM Limiting C PPM = 1-25 


32 


8.7E-3 




4.2E-2 






64 


4.5E-3 


0.9 


3.4E-2 


0.3 




128 


2.0E-3 


1.2 


2.7E-2 


0.3 




256 


8.9E-4 


1.2 


2.1E-2 


0.3 


VL+PPM Limiting C VL , C PPM = 1-25 


32 


7.7E-3 




4.1E-2 






64 


4.1E-3 


0.9 


3.3E-2 


0.3 




128 


1.8E-3 


1.2 


2.7E-2 


0.3 




256 


8.0E-4 


1.2 


2.1E-2 


0.3 
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A th Order, Square Wave 




N ceU 


L x 


Ri 


Loo 


Rod 


No Limiter 


32 


1.1E-1 




4.5E-1 






64 


7.6E-2 


0.5 


4.7E-1 


-0.1 




128 


4.5E-2 


0.8 


4.8E-1 


0.0 




256 


2.5E-2 


0.8 


4.9E-1 


0.0 


PPM Limiting C PPM = 


32 


9.1E-2 




4.2E-1 






64 


5.2E-2 


0.8 


4.3E-1 


0.0 




128 


3.0E-2 


0.8 


4.4E-1 


0.0 




256 


1.7E-2 


0.8 


4.5E-1 


0.0 


PPM Limiting C PPM = 1-25 


32 


9.0E-2 




4.2E-1 






64 


5.2E-2 


0.8 


4.3E-1 


0.0 




128 


3.0E-2 


0.8 


4.4E-1 


0.0 




256 


1.7E-2 


0.8 


4.5E-1 


0.0 


VL+PPM Limiting C VL , C PPM = 1-25 


32 


8.0E-2 




4.0E-1 






64 


4.6E-2 


0.8 


4.1E-1 


0.0 




128 


2.7E-2 


0.8 


4.2E-1 


0.0 




256 


1.6E-2 


0.8 


4.3E-1 


0.0 
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Qth Order, Gaussian Wave 






u 


Ri 


Loo 


Roo 


No Limiter 


32 


5.0E-2 




2.6E-1 






64 


1.3E-2 


2.0 


9.8E-2 


1.4 




128 


2.0E-3 


2.7 


1.8E-2 


2.4 




256 


2.6E-4 


2.9 


2.5E-3 


2.9 


PPM Limiting C PPM = 


32 


6.7E-2 




4.4E-1 






64 


2.3E-2 


1.5 


2.3E-1 


0.9 




128 


5.3E-3 


2.1 


8.7E-2 


1.4 




256 


9.2E-4 


2.5 


2.7E-2 


1.7 


PPM Limiting C PPM = 1-25 


32 


4.1E-2 




2.9E-1 






64 


1.1E-2 


1.9 


9.7E-2 


1.6 




128 


2.0E-3 


2.5 


1.8E-2 


2.4 




256 


2.6E-4 


2.9 


2.5E-3 


2.9 


VL+PPM Limiting C VL , C PPM = 1-25 


32 


4.4E-2 




3.2E-1 






64 


1.2E-2 


1.9 


1.1E-1 


1.5 




128 


2.0E-3 


2.6 


1.8E-2 


2.6 




256 


2.6E-3 


2.9 


2.6E-3 


2.9 
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Qth Order, Semi-Circle Wave 




N ce ll 


Li 


Ri 


Loo 




No Limiter 














32 


8.4E-3 




4.2E-2 






64 


3.7E-3 


1.2 


3.4E-2 


0.3 




128 


1.6E-3 


1.2 


2.7E-2 


0.3 




256 


7.4E-4 


1.1 


2.2E-2 


0.3 


PPM Limiting C PPM = 














32 


7.5E-3 




3.6E-2 






64 


3.4E-3 


1.2 


3.0E-2 


0.3 




128 


1.4E-3 


1.2 


2.4E-2 


0.3 




256 


6.3E-4 


1.2 


1.9E-2 


0.3 


PPM Limiting C PPM = 1-25 














32 


7.3E-3 




3.7E-2 






64 


3.2E-3 


1.2 


3.0E-2 


0.3 




128 


1.4E-3 


1.2 


2.4E-2 


0.3 




256 


6.1E-4 


1.2 


1.9E-2 


0.3 


VL+PPM Limiting C VL , C PPM = 1-25 














32 


7.1E-3 




3.6E-2 






64 


2.9E-3 


1.3 


2.9E-2 


0.3 




128 


1.2E-3 


1.3 


2.3E-2 


0.3 




256 


5.0E-4 


1.2 


1.9E-2 


0.3 
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Qth Order, Square Wave 




N ceU 


L x 


Ri 


Loo 


Roo 


No Limiter 


32 


9.8E-2 




4.1E-1 






64 


5.6E-2 


0.8 


4.2E-1 


0.0 




128 


3.2E-2 


0.8 


4.3E-1 


0.0 




256 


1.9E-2 


0.8 


4.4E-1 


0.0 


PPM Limiting C PPM = 


32 


7.8E-2 




3.9E-1 






64 


4.5E-2 


0.8 


4.1E-1 


0.0 




128 


2.6E-2 


0.8 


4.2E-1 


0.0 




256 


1.5E-2 


0.8 


4.3E-1 


0.0 


PPM Limiting C PPM = 1-25 


32 


7.7E-2 




3.9E-1 






64 


4.4E-2 


0.8 


4.0E-1 


-0.1 




128 


2.6E-2 


0.8 


4.2E-1 


0.0 




256 


1.5E-2 


0.8 


4.3E-1 


0.0 


VL+PPM Limiting C VL , C PPM = 1-25 


32 


6.5E-2 




3.6E-1 






64 


3.7E-2 


0.8 


3.8E-1 


-0.1 




128 


2.1E-2 


0.8 


4.0E-1 


-0.1 




256 


1.2E-2 


0.8 


4.2E-1 


-0.1 
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3.3 Comparing Conventional and Extremum-Preserving Limiters 

This section compares the conventional limiter and extremum-preserving limiter when used 
in PPM for the one-dimensional scalar advection tests. The conventional limiter uses 2 nd 
order van Leer limited differences to compute spatially 4 th order accurate extrapolated 
face-centered values [9j [TP] . The extremum-preserving limiter employs 6 th order accurate 
piecewise parabolic reconstruction without using the van Leer limited differences [ffj. As 
expected, the extremum-preserving limiters significantly reduce the error for the Gaussian 
wave (G), reduce the error only slightly for the semi-circle wave (SC), and does not reduce 
the error for the square wave (S). This comparison is shown in Figures [1] and [2j 



Conventional Limiter 



N C ell G L\ 


Ri 


G Lqo 


R<x> 


SC L x 


Ri 


S Li 


Ri 


32 7.6E-2 




4.8E-1 




7.8E-3 




8.4E-2 




64 2.7E-2 


1.5 


2.7E-1 


0.8 


4.3E-3 


0.9 


4.8E-2 


0.8 


128 7.7E-3 


1.8 


1.0E-1 


1.4 


1.9E-3 


1.2 


2.8E-2 


0.8 


256 1.3E-3 


2.6 


3.1E-2 


1.7 


8.3E-4 


1.2 


1.6E-2 


0.8 




Extremum-Preservinj 


1 Limiter (C = 1.25) 












Ncell G L\ 


Ri 


G Loo 


Roo 


SC Li 


Ri 


S Lx 


Ri 


32 4.1E-2 




2.9E-1 




7.3E-3 




7.7E-2 




64 1.1E-2 


1.9 


9.7E-2 


1.6 


3.2E-3 


1.2 


4.4E-2 


0.8 


128 2.0E-3 


2.5 


1.8E-2 


2.4 


1.4E-3 


1.2 


2.6E-2 


0.8 


256 2.6E-4 


2.9 


2.5E-3 


2.8 


6.1E-4 


1.2 


1.5E-2 


0.8 



3.4 Sensitivity to Cvl-, Cppm 

Using a Gaussian wave as a test problem, the sensitivity in the extremum-preserving limiter 
was analyzed by plotting ln(Ax) versus ln(Loo) for multiple values of Cvl = Cppm between 
1.25 and 5. The extremum-preserving limiter is insensitive to the limiting constants Cvl, 
Cppm, where Loo changed by 4% when comparing Cvl, Cppm = 1-25 and 5. 
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Figure 2: Extremum-Preserving vs Conventional Limiter. Square Wave, N ce u = 128, t — 10 
periods, CFL = 0.2, A = 1. 
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4 Algorithms for Nonlinear Systems of Hyperbolic Con- 
servation Laws 

Higher-order limiting can be extended to nonlinear systems of conservation laws. In one 
spatial dimension, these systems having the following form: 

— — h -7- — = 0, U(x, t) G R N , F = F(U), F : R N -> R N , (104) 
at ox 

where the quasilinear form is: 

A has N real eigenvalues corresponding to N linearly independent eigenvectors: 

Ar k = \ k r k , l k A = \ k l k , (106) 

r fc = r fc (f/f), l k = l k (U k ). (107) 

Like the limiting done for scalar advection, each of the variations (VL, PPM, VL+PPM) 
can be implemented for nonlinear systems of hyperbolic conservation laws. The limiters are 
either applied componentwise to the primitive variables or applied to one characteristic field 
at a time. The former approach was implemented in Chombo (adaptive mesh refinement in- 
frastructure for solving a wide variety partial differential equations) [15] and Athena (magne- 
tohydrodynamics code for astrophysical applications) [16J. Results for extremum-preserving 
limiters applied to nonlinear systems of hyperbolic conservation laws are shown in [T6J. In 
particular, [16] presents the Shu-Osher shock entropy wave interaction [17J, where the nu- 
merical solution obtained using an extremum-preserving limiter is comparable to WEN03 
and WEN05 [18] . These results are replotted below in Figures [3] and 0] with the permission 
of the authors. Figures [5] and [6] show results from the Woodward-Colella ramp problem 
[j~9] that were computed by Chombo using the unsplit PPM method of Miller-Colella [10] 
and the extremum-preserving limiters of Colella-Sekora [11]. This calculation employed a 
64 x 16 base grid with two levels of adaptive mesh refinement such that there is a factor 
of four between each level. Thus, there is an effective resolution of 1024 x 256. The two 
figures show material density for the entire domain and the double Mach region, respectively. 
From these plots, it is clear that extremum-preserving limiters are robust enough to handle 
multidimensional shocks. 
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N .. = 800 Conven 

cell 

o N „ = 200 Conven 
o cell 




-0.5 



0.5 



Figure 3: Conventional Limiter for the Shu-Osher Shock Entropy Wave Interaction. Com- 
puted using Athena. 
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Figure 4: Extremum-Preserving Limiter for the Shu-Osher Shock Entropy Wave Interaction. 
Computed using Athena. 
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Figure 5: Extremum-Preserving Limiter for the Woodward-Colella Ramp Problem. Material 
density for the entire domain. Computed using Chombo. 




Figure 6: Extremum-Preserving Limiter for the Woodward-Colella Ramp Problem. Material 
density for the double Mach region. Computed using Chombo. 
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